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Analysis of a 100 Year Geodetic Record from Northern California 

Lewis E. Gilbert^, John Beavan, and Chris Scholzt 
Lamont-Doherty Geological Observatory 
t also at Columbia University 

A geodetic network which spans the region between San Francisco and Lake 
Tahoe has been measured 5 times completely with triangulation in 1880, 1922, 

1929, 1948 and 1963. A resurvey with GPS in 1991 allows the formation of 1 
coseismic and 4 interseismic epochs. The data from this network provide a 
unique opportunity to examine the temporal and spatial evolution of the strain 
field associated with the 1906 San Francisco earthquake in particular and with the 
Pacific - North American plate boundary in general. Calculations of strain rate 
from the network data lead to following conclusions: 1) There is no resolvable 
(>0.05 p.rad/yr) strain in the between Sutter Buttes and the Sierra Nevada. 2) 
Throughout the time since the 1906 earthquake a region extending at least as far 
east as the westernmost Great Valley has been undergoing deformation related to 
Pac:Nam interaction and the associated earthquake cycle. 3) In the time and space 
of overlap, our results agree with those from the USGS trilateration data. Both 
data sets indicate that strain must be accumulating to the east of Vaca. 4) The San 
Andreas discrepancy cannot be accommodated in the Great Valley at the 1 sigma 
level of our results. It is possible to absorb it in that region at the 2 sigma level. 

5) Strain rate is elevated in the years following the earthquake and decays slowly 
with time. It is possible that the rate in the Coast Ranges increases until around 
1950 and then decays. With the exception of one epoch, strain rate in the Coast 
Ranges is consistently fault parallel, shows no sign changes and is consistent 
with monotonic strain accumulation. 

Introduction 

History of Models of the Earthquake Cycle 

In the English literature the notion of an earthquake cycle seems to have been introduced by 
Gilbert [1884] in relation to faulting in the Wasatch. Gilbert restricted his comments to the normal 
faulting of the Great Basin where large vertical offsets across normal faults are a first order geomorphic 
feature. In this setting the accumulation of large relative offsets, which is implied by repeated 
earthquakes was not conceptually difficult within the fixist world view of the time. After the great San 
Francisco earthquake of 18 April 1906, Reid [1910] extended the idea of slow accumulation and rapid 
release of elastic strain to strike-slip faults in what has come to be known as the elastic rebound 
hypothesis. Simply stated, this hypothesis asserts that earthquakes are the periodic (or at least episodic) 
release of elastic strain which has accumulated slowly during an inter-event period. A direct corollary is 
that future events will occur when the amount of strain released in the current event accumulates once 
more; thus if the strain accumulation can be monitored, earthquakes can be predicted. That prediction 
would actually be more complicated was realized by Gilbert [1909]. 

In our current historical framework the accumulated horizontal offset implicit in the elastic 
rebound hypotheses as applied to a strike-slip regime provides a mechanism for achieving large 
horizontal displacements between adjacent plates; at the time of its formulation, such offsets were an 
anathema (at least to American geologists). Large offsets could be tolerated locally, but not globally 
(see, for example, Willis [1927, p.38]). The fixist world view ol the his day led Reid to fix the relative 
positions of points at some distance from the fault (points W and E, Figure 1); thus the displacements 
close to the fault die out at some distance from the fault 
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Figure 1: Reid's diagram [Reid, 1910, Figure 12] showing the distribution of 
displacements (heavy line, WAOCE) just prior to the 1906 earthquake. The Farallon 
Islands and Mount Diablo are at the positions of A and C respectively. Reid states that 
”it is evident that the displaced area must have some limit" and represents this by having 
the displacement curve return to the zero line (WOE). 


Based on the facts that the width of the 1906 deformed zone was narrow relative to its length aiid 
that there was little evidence of extension or compression, Reid concluded that the forces responsible for 
the earthquake must have been shears concentrated close to the fault. Thus regional compressional or 
extensional forces were dismissed, as were shears acting at a distance. By the process of elimination, 
Reid postulated basal tractions as the source of the disturbing force. He asserted that those tractions 
must be the result of a northerly sub-crustal flow extending east from the fault to some distance beyond 
Mount Diablo and a southerly return flow extending to the west [Reid, 1910, p.27]. He did not attempt 
to explain the source of the requisite energy beyond vague reference to isostasy. Reid s rather awkward 
(by today's standards) mechanism is a fairly direct result of the requirement (self-evident in his day) that 
W and E (Figure 1) not move relative to each other. 

A long hiatus in thinking concerning the earthquake cycle followed Reid's work. That hiatus 
ended with the emergence of a more mobile world view in the form of plate tectonics. In the new 
paradigm, the accumulation of large offsets across a fault such as the San Andreas (S AF) is natural an 
is welcomed as a characteristic of the elastic rebound model; furthermore, points W and E in Figure 1 are 
allowed to move relative to each other, which obviates Reid's driving mechanism. 

While Reid did calculate the magnitudes of the forces necessary to produce the 1906 earthquake, 
the core of his thinking was avowedly qualitative (e.g., Reid puts no scale on Figure 1, and comments 
further along in his development that "accuracy in small details" is not insisted upon). In the current 
climate such qualitative thinking is much less in evidence and ideas concerning the earthquake cycle are 
usually presented in the form of quantitative models which attempt to approximate the physical nature ot 
the materials involved. These calculations have evolved from the elastic half-space models ot Savage 
and Hastie [1966] through layered linear elastic and viscoelastic models [Li and Rice, 1987] to non- 
linear viscoelastic models [Reches, 1992 (in press)]. 

All currently available theorizing related to the earthquake cycle (and the implicit desire to predict 
the next one) has assumed that earthquakes are roughly periodic and that the rheology important to 
earthquakes can be written down. Based on offsets measured after the 1906 San Francisco earthquake 
and displacement rates for the region in the late 1800's, Reid deduced a roughly 100 year inter-event 
time for the 1906 segment. Reid's estimate was certainly close to the mark and most modem models ot 
the earthquake cycle assume periodicity with a time scale of the order 100 years. It follows then that 
models of the earthquake cycle produce signals which develop on 100 year time scales and that 
constraint of those models will require data sets with similar temporal extent. Such data sets are rare and 
this paper presents strain calculations from the only known such data set in North America. 

In addition to constraining existing models the analysis presented below provides a valuable 
opportunity to expand our thinking concerning the earthquake cycle. It is quite likely that earthquakes 
are not periodic and in fact they may cluster [Kagan and Jackson, 1991]. Such a condition would render 
all periodic models of the earthquake cycle obsolete. Whether or not they cluster, it is clear that a 
sophisticated understanding of the cycle will require data with time scales of at least 100 years. 


The Historical Network 

Much of Reid's conceptualizing about the nature of the earthquake cycle was aided by a 
comparison of geodetic measurements made immediately after the 1906 event with those made in the 
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latter part of the 19 th century. The initial observations were made as part of the establishment of a 
transcontinental geodetic control network. This network crossed the Sierra Nevada at the latitude of 
Lake Tahoe and continued across the Great Valley to San Francisco. First order spacing resulted in the 
establishment of 3 polygons in the roughly 300 Ion between the peaks of the Sierra on the east and the 
San Andreas Fault (SAF) on the west. The Primary Arc between San Francisco and Lake Tahoe was 
established in the years between 185 1 and 1887 (Figure 2). 

After the initial observations and adjustment of the Primary Arc, San Francisco was destroyed by 
the 1906 earthquake and subsequent fire. Surface rupture of the order of 3 meters was observed dong 
the entire 280 km length of the rupture exposed on land. Geodesists of the U.S Coast and Geodetic 
survey realized that the permanent offsets associated with the earthquake would have disturbed their 
control network in the region. In order to reestablish the control network it would have to be resurveyed 
and the control points would have to be readjusted. Hayford and Baldwin [1908] report that it was 
realized that the disturbance extended beyond the easily recognizable surface rupture to a distance ot 
"many miles on each side" of the fault and that the resurvey would therefore need to cover a wide belt 
Initial triangulation to repair the network was completed in the interval between July 1906 and July 
1907 That triangulation covered a region 270 km long and 50 km wide. Most of the markers in this 
survey were of secondary and tertiary order. A least-squares fit of station position to the observed angles 
was performed and station displacements were calculated by differencing the pre- and post- earthquake 
adjusted station positions. In the course of the re-adjustment of the pre-earthquake data, evidence was 
found for a disturbance of the network during intermediate times (the disturbance is tentatively attnbuted 
to the 1868 earthquake). When all of the adjustments were finished, Hayford and Baldwin [1908 
reported that in addition the obvious displacements of 1906, another episode of movement had occurred 
around 1868. 

Hayford and Baldwin present their results as a map showing two episodes of movement; the hrst 
attributed to the 1868 earthquake and the second attributed to the 1906 event. In the light of what we 
know today some of the details of Hayford and Baldwin's [1908] map appear rather strange (e.g., Santa 
Cruz moves to the southeast in the 1906 earthquake); however a rather clear picture of the first order 
deformation associated with the SAF was understood. The importance and value of the geodetic 
network was recognized in the Lawson report [1908] and a call was made for deliberate investigation 
extending through years and decades and conducted on a wisely planned program [Lawson, 1908, p. 

The immediate post-earthquake surveys had been restricted to a region between the Farallons and 
Mount Diablo (A and C in Figure 1). It is clear that in Reid's thinking (Figure 1) displacements extend 
well beyond those points, and after the publication of Hayford and Baldwin [1908] it was reahzed that 
the new triangulation should have been started further to the east than Mount Diablo and Mocho [Bowie, 
1924], During 1922 stations from Mount Lola and Round Top to Ross Mountain and Mount Tamalpais 
were occupied. In 1923 work was begun at Loma Prieta and Mocho and continued to the south. A.t the 
end of the 1923 field season work was suspended due to the first of many funding difficulties associated 
with this network [Bowie, 1924]. At that time there was still a small gap between observations in 
southern California and those from the north. A preliminary report [Bowie, 1924] was issued based on 
the data through the 1922-23 work. In that report the magnitude of displacements increase dramatically 
to the south. When the southern and northern observations were linked it became clear that some the 
displacements were due to accumulated triangulation error and it was necessary to readjust all of the data 
west of Meades Ranch (which is -1000 km to the east in Kansas) incorporating Laplace azimuths ( ... 
[it was possible to do this] in a short time and with a relatively small amount of effort A dozen or more 
mathematicians were able to work simultaneously on the western net, and in 15 months the readjustment 
was completed..." [Bowie, 1928, p.15]). A final report was issued in 1928 [Bowie, 1928]. In that 
report displacements which were considered larger than the error were restricted to stations within about 
40 km of the 1906 rupture. 

In the solutions presented by Bowie it is necessary to hold some stations fixed and his solutions 
tend to show displacements increasing with distance from the fixed points. If the quantity under 
consideration is displacement in a region undergoing regional strain then a solution with fixed points will 
show displacement increasing (as the strain is integrated) with distance from the fixed points. Despite 
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Reid's work on the elastic strain accumulation, Bowie does not seem to have recognized this signal in 
his results. That signals might accumulate over time must have been understood as the network was 
reobserved with complete triangulations centered upon the years 1929, 1948, and 1963. 

Recent Measurements 

The elastic rebound theory enjoyed a rebirth with the emergence of plate tectonics but interest in 
the geodetic network most suited to its testing apparently died. While spotty measurements at stations 
within the network were made in the 1970's, the last complete measurements were completed in 1963 
and until our proposal to remeasure with GPS, no further measurements were planned. The lack of 
interest in the network in the early history of plate tectonics can be understood because during that time 
deformation was thought to be concentrated along narrow zones at all plate boundaries. Indeed Bowie s 
1928 report concludes that most of the measurable deformation associated with the 1906 event was 
concentrated with 25 miles of the fault. Some of his adjustments suggest motion outside of that distance 
(at Vaca and Monticello in particular) but Bowie felt that those motions were not large enough to be 
differentiated from measurement noise. The motions predicted by plate tectonics were very small 
compared to the resolving power of triangulation. 

In the early 1970's trilateration techniques with at least an order of magnitude better resolution 
over triangulation were developed by Savage and his colleagues at the USGS and line length networks 
were established throughout the San Francisco Bay area to monitor deformation associated with the 
SAF. Those networks extend to ~50 km from the fault and show definite deformation. In the North 
Bay network [Lisowski et al., 1991; Prescott and Yu, 1986] the eastern-most station is Vaca and it is 
clear from the displacement [Prescott and Yu, 1986] and velocity [Lisowski et al, 1991] profiles that 
have been presented that deformation extends beyond the eastern edge of the USGS networks. Both of 
those profiles have gradients which show no signs of approaching zero at their eastern-most stations. 

Displacement rates calculated from geodetic and geologic data can account for about 35 mm/yr of 
relative motion between North American and the Pacific. Global plate motion models [DeMets et al., 
1990; Minster and Jordan, 1984] predict relative motions between the two of the order of 50 mm/yr (the 
currently popular estimate seems to be ~48 mm/yr directed ~N35W) [DeMets et al., 1987]. The 
difference between the geodetic estimates of the plate motion and the model predictions came to be 
known as the San Andreas Discrepancy (SAD). This difference posed a serious threat when it was first 
discovered; however as plate motion models have improved and our ideas of plate boundaries have 
become more accommodating, the discrepancy has shrunk and its existence has become a source of new 
understanding rather than of difficulty [Argus and Gordon, 1991]. 

The maturation of very long baseline interferometry (VLBI) in the 1980's has contributed 
significantly to our understanding and to the resolution of the San Andreas discrepancy. Results from 
those measurements suggest that rather than a narrow boundary along the SAF, the Pac:Nam boundary 
is a mega-shear zone which includes the entire western US [Kroger et al., 1987; Ward, 1988, Ward, 
1990]. (It is interesting to note that this was first suggested by Atwater [1970] in the same paper which 
established the predictive potential of plate tectonics.) A large regional shear requires relaxation of the 
idea of narrow plate boundaries; this is complimented by allowing a large portion of the San Andreas 
discrepancy to be taken up by deformation in the Basin and Range. The most recent attack on the 
problem is presented by Argus and Gordon [1991]. They break Pac:Nam relative motion into three 
components: 1) SAF motion, 2) relative motion between a rigid block consisting of the Sierra Nevada 
and the Great Valley (Sierra Nevada block) and North America and 3) the remaining discrepancy. 
Motion of the Sierra Nevada block is basically the integral of Great Basin deformation. VLBI data are 
inverted to constrain these motions and the resulting discrepancy is reduced to 8 mm/yr directed N15E. 

Summary 

All of the issues touched on above suggest that a data set of regional extent and long history 
would contribute greatly to our understanding of the earthquake cycle. Data from the northern California 
Primary Arc installed by the U.S. Coast and Geodetic survey in the late 1800's form such a data set. 
We have extended that network's history beyond the 100 year mark by repairing and reobserving it 
using modem techniques (GPS). Below we present the results of our new measurements and a 
reanalysis of the historical data. As the initial measurements in the data presented here provided Reid 
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with information fundamental to his formulation of the elastic rebound theory, so too the data presented 
below may again inspire a unique insight to the workings of Earth. 

Data Analysis 


Subdivision of the data 

In the current study we have restricted ourselves to first order observations made from, and to, 
the first order control points east of the SAF (Figure 2). This allows us to use only the best possible 
data and to concentrate on shear strain accumulation without introducing an explicit scale mto the 
problem. A major ramification of this decision is that our first post- 1906 measurements are from 1922. 
We are currently integrating the lower order near-fault data (collected in 1906-1907) into our analysis 
and the results of that work will be presented in a future publication. 

In calculating strain it is necessary to assume that the strain rate is constant over the time and area 
of consideration. Our data set includes a major earthquake and spans a region which includes three 
geographic provinces. It is not reasonable to assume that strain is uniform over large portions of the set 
and we have subdivided the data into sets for which a uniform strain assumption is more reasonable. 

The data have been subdivided geographically as indicated in Figure 2. At the finest scale 
possible the data break naturally into three geographic provinces: 1) the Coast Ranges which is bordered 
by the stations Ross Mountain, Mount Tamalpais 2, Mount Diablo, Mount Vaca, Monticello, and Mount 
Helena- 2) the Great Valley which is bordered by the stations Mount Diablo, Mount Vaca, Monticello, 
Mount Helena, Marysville Buttes, and Pine Hill; and 3) the Foothills which is bounded by Marysville 
Buttes, Pine Hill, Mount Lola and Round Top. 

It is possible to break the data set into 5 temporal subsets; a coseismic epoch which includes the 
1906 event and 4 interseismic epochs. The coseismic epoch is formed from all of the pre-earthquake 
data (observations between 1878 and 1903) and the first complete survey following the earthquake 
(1922). It commonly took several field seasons to observe the entire network with triangulation; thus 
interseismic epochs are formed from complete triangulation surveys centered roughly on the years 1922, 
1929, 1948, 1963 and from the GPS data collected in 1991. 

Processing of the triangulation 

In each of the regions and for each of the epochs, strain rate has been calculated using the 
simultaneous reduction technique of Bibby [1982] as coded in the National Geodetic Survey program 
Dynap [ Snay and Drew, 1988]. If the assumptions of temporally and spatially uniform strain are not 
met by the data in question, the calculated strain represents some average over the varying part and the 
associated error will be larger to the extent that homogeneity is not a good approximation. 

In all except our last (GPS) epoch our observations consist almost exclusively of triangulation 
measurements (we have a limited number of astronomical azimuths). In order for the calculations to be 
non-singular it is necessary to provide some information on the scale and orientation of the network. 
For each strain calculation the azimuth and distance between two stations at the beginning and end times 
are included as observations. In cases where astronomic azimuths are available, those are used; 
otherwise an azimuth is calculated using the a priori geodetic coordinates. In all cases it is necessary to 
calculate a distance based upon the a priori coordinates. The end point stations are chosen in such a way 
that the constrained line is roughly parallel to and as far as possible from the SAF. In addition to the 
scale and orientation constraints, we must fix the coordinates of one station. The details of these 
constraints have little effect on the final results. 

As noted above, our observations are limited to horizontal directions; thus the only strain 
components that we can meaningfully compute are horizontal shears. Conventionally these are presented 
as the "engineering" shear strain components Yi and 75 : 

Yi =£ 22 -£11 
72 = £12 + £21 

where Ea are components of the calculated strain tensor and the 1 and 2 directions are North and East 
respectively. In this study we have rotated these components 15 degrees to the east so that 71 represents 
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right lateral shear strain across vertical planes oriented N30W (roughly parallel to the SAF) and Y2 
represents left lateral shear strain across planes oriented N15E (Figure 3). 



Figure 3: The shear strain components in a fault parallel and fault perpendicular 
reference frame. Yi represents right lateral shear strain across vertical planes oriented 
N30W (roughly parallel to the SAF) and 72 represents left lateral shear strain across 
planes oriented N15E. 


Processing of the GPS . f 

The GPS data were processed using the Bernese version 3.3 software. With the exception ot 
day 080 when no Marysville data was collected, all NCPAR data was differenced with Marys vie 
Buttes. On day 080 data were differenced with Mount Diablo. Cycle slips were removed from the 

single difference files. „ , . . 

In ordere to improve the orbits, CIGNET data [Chin, 1988] from Kokee Westford, Richmond, 
and Mojave were also processed. Those files were differenced with Marysville Buttes and with varying 
degrees of success, cycle slips were removed from the single differences. Removing cycle slips from 
the lines with Minimac receivers (Westford, Richmond, and Mojave) was very time consuming 
eventually a point was reached where the trade off between further effort and improved quality did not 
warrant more effort [Larson andAgnew, 1991a]. 

After cycle slips were removed from the CIGNET data, daily ionosphere solutions for t e 
nosition of Marysville Buttes and for the satellite orbits were calculated with the CIGNET station 
positions fixed. The resulting positions for Marysville Buttes were averaged to arrive at a mean position 
for that station. The scatter of the daily positions with respect to this mean are presented in Figure 4 (a). 
It is clear that the scatter in the north component is smaller than that in both the east and up components 
r Larson and Asnew, 1991b]. Histograms of the scatter are presented m Figure 4 (b). Given the small 
number of data (9) the scatter is acceptably symmetric about zero. If we take 4000 km as the average 
baseline in this solution and 10 cm as the scatter in the repeats, this represents a precision of 2.5 parts in 
108 . The Mojave baseline is only 600 km and yields a precision somewhat worse. 
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Figure 4: Scatter of Marysville Buttes position when calculated using fixed CIGNET 
stations and free orbits. In (a) both North and Up scatters are plotted against East. In (b) 
histograms of the scatter are presented; tick increment on the vertical axis is 1. 

The calculated mean position for Marysville Buttes was then held fixed along with the CIGNET 
stations and daily ionosphere free solutions for the NCPAR station positions and orbits were calculated. 
To investigate the quality of the daily solutions, the calculated interstation vectors were adjusted using 
Dynap. On the basis of the adjustment several of the baselines measurements were discarded and the 
remaining vectors were adjusted again. The quality of the remaining solutions is illustrated in Figure 5. 
In that plot we show the standard deviation of the daily GPS baseline lengths about the adjusted length. 
The quality of the solutions is fairly independent ol baseline length and the standard deviations have an 
rms scatter of about 1.3 cm. 
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Also shown in Figure 5 is a measure of precision (the standard deviation divided by the baseline 
length) plotted against baseline length. In the triangulation the precision is about 1CP. This is far less 
precise than our daily solutions. For the purposes of comparing the GPS with the triangulation it we 
decided against attempting to fix the ambiguities in the GPS solution. 

To calculate the strain in our final epoch, the daily solutions were used to calculate the dx,dy and 
dz components of the measured baselines. Based on the repeatabilities (Figure 5) a conservative a priori 
sigma of 2 cm was assigned to each component. Finally the GPS data and the triangulation data were 
combined using Dynap to calculate strain in epochs ending in 1991. 


Repeatability of Daily Solutions 9 Standard Deviation 

within the NCPAR network O Std. Dev./Baseline Length 



0 50 100 150 200 250 


Baseline Length (km) 

Figure 5: Repeatability of the daily solutions. These ionosphere free solutions are 
calculated by holding the positions of the CIGNET stations and Marysville Buttes fixed 
and solving for the orbits and NCPAR station coordinates. Standard deviations are 
calculated relative to the adjusted baseline length. 


The results of this study are presented in Figure 6 and tabulated in Table 1. In Figure 6 the 
results for each region (Figure 2) and each strain component are plotted against time. The time for each 
point is the mid-point of the epoch for that calculation. Distance from the fault increases from left to 
right in the figure with plots from the Coast Range leftmost, Great Valley center and Foothills on the 
right. From top to bottom the components yi, 72, maximum shear strain rate and azimuth ot minimum 
extension are presented. Error bars represent the calculated 1 sigma error. In addition to representing 
scatter in the data those error bars will also reflect the extent to which the region fails to be adequately 
described by a uniform strain model; thus if strain rate varies significantly over the time or area of our 
calculations, the calculated error will be larger than that expected simply from observational error. 

In considering Figure 6 two things to keep in mind are: 1) Yi and 72 have been rotated by 15 
degrees to the east so that 71 is roughly parallel to the SAF (Figure 3). Signal in 72 represents strain 
accumulation which is oblique to the plate boundary. 2) Results for the first epoch include data from 
-1880 to 1922. That epoch includes pre- co- and post-seismic signal. Interpretation of that result is 
difficult and will not be attempted here. 

A schematic presentation of the results is given in Figure 7. In that figure the length of the lines 
represent the magnitude of the maximum shear rate and the orientation that of the minimum extension 
direction relative to the SAF (lines oriented at 45 degrees will produce simple shear across the SAF). 
The ellipses represent the error; both axes have length representing 1 sigma. The axis parallel to the line 
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is the error in the rate; the perpendicular axis is the error in the direction. Error ellipses overlap when the 
shear rate is within 1 sigma of zero. 

Notice that as expected within the plate tectonic paradigm and from the results of Argus [1991], 
the signal in the Foothills province always is indistinguishable from zero. Strain rate there vanes 
smoothly (decreasing from a maximum immediately after the earthquake), but the errors associated with 
the calculated results are quite large relative to the signal; thus the first result of this study is that strain 
rate in the Foothills region is smaller than we can detect with any confidence. An upper limit on its 
current magnitude is about 0.05 ^rad/yr (Table 1). For the remainder of this paper we will consider 
strain in the Foothills region to be beyond the limit of detectable effects associated with the plate 
boundary (strain accumulation related to Pac:Nam relative motion and strain migration associated with 
the 1906 earthquake). 

The results from the Coast Ranges in the epoch 1948-1963 do not fit into any monotonic scheme 
of an earthquake cycle. The data from this epoch contain no obvious blunders or extreme outliers and it 
appears that they can not be simply discarded. If the data from 1948-1963 and 1963-1991 are combined 
(open dots in Figure 6) the results are consistent with monotonic strain accumulation and more in line 
with the earlier results; thus the anomaly is restricted to the combination of the 1948 and 1963 data. In 
the discussion that follows we will address the combined result (1948-1991) and we will return to 
speculate on the 1948-1963 result at the end of the paper. 

With the above caveats in mind (no detectable strain in the Foothills and ignore the Coast Range 
1948-1963 result), several first order observations concerning the evolution of the strain field related to 
the 1906 earthquake and the Pac:Nam boundary can be made: 1) Strain associated with the 1906 
earthquake is not restricted to the Coast Ranges. Maximum shear rate in the Great Valley is high after 
the earthquake and decreases monotonically through to the present. Inspection of Figure 7 demonstrates 
graphically that the results between 1922 and 1948 for both the Coast Ranges and Foothills are 
consistent with strain primarily parallel to the SAF. 2) Throughout the Coast Ranges, fault parallel strain 
decays over a period of at least 40 years following the earthquake. During this time in this region there 
is no significant strain oblique to the plate boundary. 3) The azimuth results from the Great Valley 
suggest that the orientation of the strain field in that region rotates with time. 


Coast Ranges 


Years 

Yi 

72 

1880- 1922 

-0.12 ± 0.08 

-0.21 ± 

1922 - 1929 

0.52 ± 0.74 

-0.17 ± 

1929 - 1948 

1.00 ± 0.42 

0.06 ± 

1948 - 1963 

-0.25 ± 0.37 

-1.34 ± 

1963 - 1991 

0.41 ± 0.09 

-0.02 ± 

1922 - 1948 

0.53 ± 0.13 

-0.02 ± 

1948- 1991 

0.33 ± 0.05 

-0.03 ± 

1922 - 1963 

0.30 ± 0.10 

-0.02 ± 

1922 - 1991 

0.35 ± 0.03 

-0.04 ± 



Max. Shear 

0.10 

0.24 ± 0.11 

0.96 

0.55 ± 0.63 

0.42 

1.00 ± 0.43 

0.75 

1.37 ± 0.75 

0.16 

0.41 ± 0.09 

0.14 

0.53 ± 0.12 

0.08 

0.33 ± 0.05 

0.13 

0.30 ± 0.09 

0.04 

0.35 ± 0.03 


Azimuth 


74.85 

± 

8.62 

24.16 

± 

53.69 

13.31 

± 

11.83 

65.36 

± 

7.82 

16.33 

± 

11.42 

16.14 

± 

7.70 

17.83 

± 

7.20 

16.53 

+ 

11.82 

18.18 

± 

3.40 
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Great Valley 


Years 

Y * 1 

Yi 

Y2 

Max. Shear 

Azimuth 

1880-1922 

0.01 ± 0.09 

-0.22 ± 0.12 

0.22 ± 0.12 

58.95 ± 11.97 

1922 - 1929 

0.65 ± 0.52 

0.36 ± 0.62 

0.75 ± 0.54 

0.58 ± 23.21 

1929 - 1948 

0.42 ± 0.30 

-0.62 ± 0.30 

0.75 ± 0.34 

42.89 ± 9.66 

1948 - 1963 

-0.23 ± 0.28 

-0.30 ± 0.28 

0.38 ± 0.23 

78.37 ± 24.3 

1963- 1991 

-0.05 ± 0.07 

0.00 ± 0.07 

0.05 ± 0.07 

-75.61 ± 40.55 

1922 - 1948 

-0.11 ± 0.34 

-0.09 ± 0.23 

0.15 ± 0.20 

85.02 ± 70.12 

1948 - 1991 

-0.05 ± 0.04 

-0.01 ± 0.04 

0.05 ± 0.04 

-79.01 ± 23.15 

1922 - 1963 

0.01 ± 0.23 

-0.16 ± 0.17 

0.16 ± 0.18 

58.14 ± 40.18 

1922- 1991 

0.05 ± 0.03 

-0.09 ± 0.03 

0.11 ± 0.03 

46.95 ± 6.85 


Foothills 


Years 

Yl 

Y2 

Max. 

Shear 

Azimuth 

1880 - 1922 

-0.18 

± 

0.86 

-0.05 

± 

0.75 

0.19 

± 

0.69 

-83.15 

± 

147.2 

1922 - 1929 

-0.35 

± 

2.19 

0.50 

± 

1.31 

0.61 

± 

2.13 

-47.71 

± 

71.02 

1929 - 1948 

-0.08 

± 

0.37 

0.21 

± 

0.48 

0.23 

± 

0.46 

-40.08 

± 

50.99 

1948 - 1963 

0.06 

± 

0.44 

0.04 

+ 

0.46 

0.08 

± 

0.43 

-1.71 

± 

172.4 

1963 - 1991 

-0.12 

± 

0.09 

-0.04 

± 

0.13 

0.12 

± 

0.11 

-83.90 

± 

28.64 

1922 - 1948 

-0.75 

± 

1.71 

0.37 

± 

0.34 

0.83 

± 

1.66 

-61.92 

± 

20.81 

1948 - 1991 

-0.07 

± 

0.07 

-0.01 

± 

0.07 

0.07 

± 

0.07 

-79.17 

± 

30.70 

1922 - 1963 

-0.05 

± 

0.04 

-0.02 

± 

0.04 

0.06 

± 

0.05 

-87.43 

± 

20.77 

1922- 1991 

-0.05 

± 

0.04 

-0.02 

± 

0.04 

0.05 

± 

0.05 

-87.91 

± 

21.13 

Table 1: Tabulated 

1 results of 

this study. 

yi is parallel to the SAF. 

. Strain 

rates are in 

pradian/yr and the azimuth is 

that of the minimum extension direction, brrors are i 


sigma. 


Discussion 

Spatial and Temporal Breadth of the Strain Field 

VLBI results suggest that the deformation associated with Pac:Nam relative motion is spread 
across the entire western US (e.g. [Ward, 1990]); however the spatial distribution of VLBI stations does 
not allow distinction of the details. At the finer end of the spatial spectrum, the USGS trilateration 
networks provide fairly detailed information concerning the recent spatial distribution of strain within the 
Coast Ranges but those networks do not extend east of Vaca. With roughly 100 km baselines spanning 
the region between San Francisco and Lake Tahoe, our network fills a significant gap in the available 
spatial bandwidth; with an age of over 100 years our network also extends considerably our temporal 
bandwidth. 

Our most recent (1963-1991) Coast Ranges result overlaps with the USGS trilateration in both 
space and time. It is reassuring that results from the two are comparable (Table 2). As noted above, the 
trilateration results require that strain currently continues into the Great Valley east of Vaca. Our Great 
Valley results for 1963-1991 are very small yet they are certainly consistent with strain inside the region. 
The Great Valley results are averaged over a very large area; this combined with the errors m the 
triangulation render our most recent result close to our measurement threshold. Our Great Valley strain 
results could be accounted for by holding Marysville Buttes and Pine Hill fixed and moving only the 
markers at the western boundary of the zone (Vaca sits atop the first ridge coming west out of the Great 
Valley and Monticello is only a couple of ridges to the west of the great flatness of the valley). This 
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would allow the Great Valley strain to be concentrated in the westernmost portion of the region. In the 
earlier epochs strain rate in the Great Valley was much higher than it is currently; thus it appears that the 
Great Valley has been affected by strain related to the Pac:Nam boundary throughout the history ot the 

Primary Arc network. 


Lisowski et al. [1991] 
This Study 


Max. Shear Rate 
(engineering) 


Azimuth Minimum. 
Extension 


N10.6E ± 1.5 
N16.3E ± 11.4 


0.34 ± 0.02 

, 0.41 ±0.09 . _ , 

Table 2: Comparison of our 1963-1991 Coast Ranges results with those from the North 
Bay trilateration networks. The results are the same within the 1 sigma errors. 


In at least the 40 years following the earthquake the strain rate in both the Great Valley and the 
Coast Ranges was elevated relative to its current level (considerably so in the Great Valley). These rates 
appear to be highest close to the earthquake and to decay slowly (m the Coast Ranges it is possible tha 
the rate increases into the 1929-1948 epoch and does not begin to decay until after that time) (Figures R 
and 7). The character of the strain rate decay in our regional results contrasts strongly with the rapid 
decay suggested by the near-fault data compiled by Thatcher [1983]. With the available data it is not 
possible to determine whether the decay in the Great Valley is due simply to a decrease m rate within the 
region currently undergoing strain or to a contraction of the area undergoing strain. 

Regardless of the details, it should be recognized that strain related to the earthquake cycle of the 
1906 segment extends at least as far east as the westernmost Great Valley. The scenario ofleast 
astonishment seems to be that strain extends to the east of Vaca (as required by the trilateration) but is 
concentrated in the westernmost Great Valley (leading to a diminished signal in the spatially averaged 
triangulation); the strain rate is elevated through out the Coast Ranges and western Great Valley shortly 
after the earthquake and decays over a period of at least 50 years following the earthquake. 

Attempts to find the San Andreas discrepancy still fall about 8 mm/yr shon [Argus (md Gordon, 
1991] and it is possible that the remaining strain is accommodated in the Great Valley. If all ot the 8 
mm/yr is taken up as SAF parallel shear strain in the 100 km between the Coast Ranges and the Foothills 
we Jould expect strain rates of 0.08 prad/yr. We measure -0.05 ± 0.07 prad/yr of fault paralle shear in 
the Great Valley (Table 1: 1963-1991 Great Valley yi). This value is indicates strain release rather than 
accumulation and we find that our results will not accommodate the discrepancy as SAF parallel shear 
strain in the Great Valley at the 1 sigma level. It could be accommodated in that region at the 2 sigma 

level. 

Implications of the Strain Field Evolution 

Data relating to the temporal and spatial extent of the strain field are interesting to those who 
study the earthquake cycle because they provide some definition of our conceptual realm. Our hope is 
that an acceptable analog to the earthquake cycle can be found and that that analog wdl allow us to predict 
future events. Currently our conceptualizing is dominated by quantitative physical models which take 
the form of systems of equations describing the rheology of a region invovled m the earthquake cycle. 
These systems are subjected to periodic disturbances which are analogous to earthquakes and the 
evolution of the equations is studied. If an analog is to be judged a suitable, it must reproduce the 
monitored behavior of Earth during the earthquake cycle. . , 

A common characteristic of many models of the earthquake cycle is that they have a horizontal y 
layered rheology. An elastic surface layer is usually underlain by some sort ofweaker layer; in the 
model of Li and Rice [Li and Rice, 1987] this layer has a viscoelastic rheology. The Li and Rice model 
has 5 parameters; plate velocity, depth of the locked zone, and cycle period pre fixed on basis of 
geologic and seismologic evidence. The two remaining parameters, thickness ot the elastic layer and 
relaxation time of the underlying viscoelastic layer are available for adjustment in light ot surface 
geodetic data. 

In the current study we have analyzed observations over temporal and spatial scales which 
exceed those previously available by a factor of at least 2 or 3. In so doing we have introduced new 
characteristics which must be produced by acceptable models of the earthquake cycle. While our data fill 


NASA NAG 5 -799 Final Report 


- 11 - 


March 1992 


an important gap in the spatial and temporal record, they are sparse compared to the continuous profiles 
which can be produced by the models. Our data are in the form of spatial and temporal averages of 
strain rate; therefore to make a detailed comparison between our data and a model it will be necessary to 
form averages of the model predictions. A detailed comparison (forming the necessary averages) is 
beyond the scope of this paper but we will qualitatively explore the implications of the current results as 
they relate to the model of Li and Rice [Li and Rice, 1987]. 

From the results presented above (Figure 6 and 7) it is clear that, at least early in the earthquake 
cycle, the data require elevated strain rates at distances beyond 80 km from the fault. In the model of Li 
and Rice the width of the strain field is most easily related to the thickness of the lithosphere and they 
prefer values of 20-30 km for this parameter. If 20 km is chosen, strain rates beyond 70 km are always 
very low; our results do not allow this. If 30 km is chosen, elevated strain rates can be achieved beyond 
75 km very early in the cycle; thus our results tend to favor a thicker elastic layer in the scenario of Li 
and Rice. 

Our data also require that strain rates remain high for a significant portion of the earthquake 
cycle. For a given lithospheric thickness, strain rate late in the cycle can be sustained by increasing the 
viscoelastic relaxation time; thus our results tend to support longer relaxation times. It must be noted that 
increasing the thickness to broaden the field will also decrease the strain rate late in the cycle, that is, 
increasing relaxation time and thickening the elastic layer balance each other. It is not clear that the 
breadth in both time and space implied by our results can be matched by varying, within reasonable 
limits, the thickness and relaxation time parameters of the Li and Rice model. 

A possible solution to this dilemma is hinted at by Li and Rice [1987]. They note that several 
lines of evidence point to a shallowing of the locked zone during the evolution of the earthquake cycle, 
(e.g. [Tse and Rice, 1986]). They do not model this phenomenon but they note such an effect would 
tend to elevate late cycle rates and might counteract the decreasing rates present in their model. 

Speculations re: the 1948-1963 epoch 

In the discussion up to this point we have ignored the 1948-1963 Coast Ranges results and in its 
place have considered the ensemble result from 1948 to 1991. We chose to proceed in this fashion 
because the 1948-1963 Coast Ranges result is of a strikingly different character than the results which 
surround it. With the elimination of the Coast Ranges 1948-1963 point, our strain rate results from that 
region are consistently fault parallel, show no sign changes and are consistent with monotonic strain 
accumulation. This is a tidy situation, but it must be disturbed because the 1948-1963 result can not 
simply be discarded. We now return to consider it further. 

In the Coast Ranges, yi increases for 40 years after 1906 (consistent with strain accumulation) 
and then abruptly changes sign to become slightly negative; concurrent with the yi sign change, Y2 
becomes strongly negative after having fluctuating about zero since the earthquake (Figure 6). In terms 
of maximum rate and orientation, the Coast Ranges azimuth is consistent with fault-parallel right lateral 
shear strain accumulation during the first two epochs and the magnitude of the maximum shear strain rate 
increases steadily. In the 1948-1963 epoch the rate continues to increase but the orientation of the strain 
field changes dramatically from fault parallel to strongly oblique (Figures 6 and 7). The results from the 
1948-1963 Coast Ranges epoch indicate that the region was dominated by SAF oblique strain release 
during that period. 

The 1948-1963 results are germane to an important incompatibility between our data and the Li 
and Rice model. At distances greater than 1 lithospheric thickness the model predicts negative (release) 
strain rates for roughly the first 25 years (—15% of 160 year cycle). If we exclude the 1948-1963 results 
from consideration, our Coast Ranges results (those spatially comparable to the model) do not show any 
sign change at all. It is possible that the spatial and temporal averaging inherent in our results has 
masked the effect; however if this was the case we would expect our early results to be lower than they 
are. With the reintroduction of the 1948-1963 results our results include a sign change; none-the-less 
there are still some incompatibilities between our results and the model. The largest incompatibility has 
to do with the timing of the sign change. Our strain release result follows a 40 year period of strain 
accumulation; in the model strain release follows the earthquake immediately. If, for the sake of 
argument, we ignore this detail and match the 1948-1963 results with the model strain release, it is 
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